Chapter 6 Beta diversity

load("data/data.Rdata")
beta_q0n <- genome_counts_filt %>%
  column_to_rownames(., "genome") %>%
  filter(rowSums(. != 0, na.rm = TRUE) > 0) %>%
  dplyr::select_if(~!all(. == 0)) %>%
  hillpair(., q = 0)

beta_q1n <- genome_counts_filt %>%
  column_to_rownames(., "genome") %>%
  filter(rowSums(. != 0, na.rm = TRUE) > 0) %>%
  dplyr::select_if(~!all(. == 0)) %>%
  hillpair(., q = 1)

beta_q1p <- genome_counts_filt %>%
  column_to_rownames(., "genome") %>%
  filter(rowSums(. != 0, na.rm = TRUE) > 0) %>%
  dplyr::select_if(~!all(. == 0)) %>%
  hillpair(., q = 1, tree = genome_tree)

beta_q1f <- genome_counts_filt %>%
  column_to_rownames(., "genome") %>%
  filter(rowSums(. != 0, na.rm = TRUE) > 0) %>%
  dplyr::select_if(~!all(. == 0)) %>%
  hillpair(., q = 1, dist = dist)

6.1 location

6.1.1 Richness diversity

betadisper(beta_q0n$S, sample_metadata$location) %>% permutest(., pairwise=TRUE) 

Permutation test for homogeneity of multivariate dispersions
Permutation: free
Number of permutations: 999

Response: Distances
          Df  Sum Sq  Mean Sq      F N.Perm Pr(>F)    
Groups     5 0.51463 0.102926 11.076    999  0.001 ***
Residuals 86 0.79917 0.009293                         
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Pairwise comparisons:
(Observed p-value below diagonal, permuted p-value above diagonal)
               Aruba     Brazil  CaboVerde      Spain    Denmark Malaysia
Aruba                1.0000e-03 2.0000e-03 1.0000e-03 1.0000e-03    0.015
Brazil    8.1559e-05            1.3700e-01 6.8000e-02 9.0600e-01    0.151
CaboVerde 6.2922e-04 1.1254e-01            1.0000e-03 3.0000e-02    0.937
Spain     1.7363e-12 7.3505e-02 3.5940e-05            2.1000e-02    0.001
Denmark   1.3206e-07 9.0556e-01 2.9720e-02 2.3526e-02               0.051
Malaysia  9.6128e-03 1.3730e-01 9.2726e-01 4.6570e-04 5.9273e-02         
adonis2(beta_q0n$S ~ location, 
        data = sample_metadata %>% arrange(match(sample,labels(beta_q0n$S))), 
        permutations = 999) %>%
        broom::tidy() %>%
        tt()
tinytable_xckn08xbszxb3z43bsl8
term df SumOfSqs R2 statistic p.value
location 5 6.57607 0.2442751 5.559605 0.001
Residual 86 20.34468 0.7557249 NA NA
Total 91 26.92075 1.0000000 NA NA
pairwise.adonis(beta_q0n$S, sample_metadata$location, perm = 999)
                   pairs Df SumsOfSqs  F.Model         R2 p.value p.adjusted sig
1        Aruba vs Brazil  1 2.0520609 7.280171 0.19528265   0.001      0.015   .
2     Aruba vs CaboVerde  1 1.0026264 3.213494 0.09975614   0.002      0.030   .
3       Aruba vs Denmark  1 1.9513181 7.018855 0.19486614   0.001      0.015   .
4      Aruba vs Malaysia  1 1.1650328 3.678697 0.11257172   0.001      0.015   .
5         Aruba vs Spain  1 2.0433431 8.167535 0.21974917   0.001      0.015   .
6    Brazil vs CaboVerde  1 2.1090004 9.164611 0.24013374   0.001      0.015   .
7      Brazil vs Denmark  1 0.5703370 2.907951 0.09113563   0.001      0.015   .
8     Brazil vs Malaysia  1 0.9283278 3.953419 0.11996993   0.001      0.015   .
9        Brazil vs Spain  1 0.6343957 3.769479 0.11503016   0.001      0.015   .
10  CaboVerde vs Denmark  1 1.8099512 8.070070 0.22373314   0.001      0.015   .
11 CaboVerde vs Malaysia  1 1.2143821 4.593887 0.14094322   0.001      0.015   .
12    CaboVerde vs Spain  1 1.8141430 9.281722 0.24896172   0.001      0.015   .
13   Denmark vs Malaysia  1 1.0124763 4.418610 0.13629855   0.001      0.015   .
14      Denmark vs Spain  1 0.5305346 3.310771 0.10573905   0.001      0.015   .
15     Malaysia vs Spain  1 0.8450695 4.218747 0.13094076   0.001      0.015   .
#pdf("figures/beta_q0_loca.pdf",width=9, height=5)
beta_q0n$S %>%
  vegan::metaMDS(., trymax = 500, k = 2, trace=0) %>%
  vegan::scores() %>%
  as_tibble(., rownames = "sample") %>%
  left_join(sample_metadata, by = "sample") %>%
  group_by(location) %>%
  mutate(x_cen = mean(NMDS1, na.rm = TRUE)) %>%
  mutate(y_cen = mean(NMDS2, na.rm = TRUE)) %>%
  ungroup() %>%
  ggplot(aes(x = NMDS1, y = NMDS2, color = location, fill = location)) +
    geom_point(size = 4) +
    #   stat_ellipse(aes(color = beta_q1n_nmds$Groups))+
    geom_segment(aes(x = x_cen, y = y_cen, xend = NMDS1, yend = NMDS2), alpha = 0.9, show.legend = FALSE) +
    scale_color_manual(values = location_colors)+
    theme_classic() +
    theme(
      axis.text.x = element_text(size = 12),
      axis.text.y = element_text(size = 12),
      axis.title = element_text(size = 20, face = "bold"),
      axis.text = element_text(face = "bold", size = 18),
      panel.background = element_blank(),
      axis.line = element_line(size = 0.5, linetype = "solid", colour = "black"),
      legend.text = element_text(size = 16),
      legend.title = element_text(size = 18),
      legend.position = "right", legend.box = "vertical"
    ) 

#dev.off()

6.1.2 Neutral diversity

betadisper(beta_q1n$S, sample_metadata$location) %>% permutest(., pairwise=TRUE) 

Permutation test for homogeneity of multivariate dispersions
Permutation: free
Number of permutations: 999

Response: Distances
          Df Sum Sq  Mean Sq      F N.Perm Pr(>F)   
Groups     5 0.2728 0.054560 3.7079    999  0.006 **
Residuals 86 1.2654 0.014715                        
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Pairwise comparisons:
(Observed p-value below diagonal, permuted p-value above diagonal)
               Aruba     Brazil  CaboVerde      Spain    Denmark Malaysia
Aruba                1.1000e-02 1.0000e-03 1.0000e-03 1.0000e-03    0.002
Brazil    9.2383e-03            7.7300e-01 2.4400e-01 6.4900e-01    0.881
CaboVerde 7.7521e-04 7.7033e-01            3.1300e-01 8.8200e-01    0.890
Spain     1.6646e-05 2.2327e-01 3.0424e-01            3.7200e-01    0.258
Denmark   3.3254e-04 6.6322e-01 8.7642e-01 3.7291e-01               0.785
Malaysia  1.1422e-03 8.5989e-01 8.9763e-01 2.4754e-01 7.7435e-01         
adonis2(beta_q1n$S ~ location, 
        data = sample_metadata %>% arrange(match(sample,labels(beta_q1n$S))), 
        permutations = 999) %>%
        broom::tidy() %>%
        tt()
tinytable_t3832fvl61scxg31xvoo
term df SumOfSqs R2 statistic p.value
location 5 5.709562 0.221304 4.88821 0.001
Residual 86 20.090068 0.778696 NA NA
Total 91 25.799630 1.000000 NA NA
pairwise.adonis(beta_q1n$S, sample_metadata$location, perm = 999)
                   pairs Df SumsOfSqs  F.Model         R2 p.value p.adjusted sig
1        Aruba vs Brazil  1 1.3739907 4.789154 0.13766227   0.001      0.015   .
2     Aruba vs CaboVerde  1 1.0636589 3.806728 0.11603497   0.001      0.015   .
3       Aruba vs Denmark  1 1.6037782 5.813780 0.16699651   0.001      0.015   .
4      Aruba vs Malaysia  1 1.2584980 4.466640 0.13346545   0.001      0.015   .
5         Aruba vs Spain  1 1.4382789 5.558030 0.16083181   0.001      0.015   .
6    Brazil vs CaboVerde  1 1.5901673 7.027154 0.19505161   0.001      0.015   .
7      Brazil vs Denmark  1 0.6954159 3.122217 0.09719804   0.001      0.015   .
8     Brazil vs Malaysia  1 0.7313843 3.199010 0.09935119   0.001      0.015   .
9        Brazil vs Spain  1 0.4007935 1.948927 0.06297236   0.034      0.510    
10  CaboVerde vs Denmark  1 1.8198218 8.556140 0.23405480   0.001      0.015   .
11 CaboVerde vs Malaysia  1 1.3171418 6.019861 0.17695136   0.001      0.015   .
12    CaboVerde vs Spain  1 1.5415454 7.905417 0.22017337   0.001      0.015   .
13   Denmark vs Malaysia  1 1.0691578 4.970172 0.15074753   0.001      0.015   .
14      Denmark vs Spain  1 0.5766622 3.014217 0.09718823   0.004      0.060    
15     Malaysia vs Spain  1 0.6446972 3.265586 0.10444666   0.001      0.015   .
#pdf("figures/beta_q1n_loca.pdf",width=9, height=5)
beta_q1n$S %>%
  vegan::metaMDS(., trymax = 500, k = 2, trace=0) %>%
  vegan::scores() %>%
  as_tibble(., rownames = "sample") %>%
  left_join(sample_metadata, by = "sample") %>%
  group_by(location) %>%
  mutate(x_cen = mean(NMDS1, na.rm = TRUE)) %>%
  mutate(y_cen = mean(NMDS2, na.rm = TRUE)) %>%
  ungroup() %>%
  ggplot(aes(x = NMDS1, y = NMDS2, color = location, fill = location)) +
    geom_point(size = 4) +
    #   stat_ellipse(aes(color = beta_q1n_nmds$Groups))+
    geom_segment(aes(x = x_cen, y = y_cen, xend = NMDS1, yend = NMDS2), alpha = 0.9, show.legend = FALSE) +
    scale_color_manual(values = location_colors)+
    theme_classic() +
    theme(
      axis.text.x = element_text(size = 12),
      axis.text.y = element_text(size = 12),
      axis.title = element_text(size = 20, face = "bold"),
      axis.text = element_text(face = "bold", size = 18),
      panel.background = element_blank(),
      axis.line = element_line(size = 0.5, linetype = "solid", colour = "black"),
      legend.text = element_text(size = 16),
      legend.title = element_text(size = 18),
      legend.position = "right", legend.box = "vertical"
    ) 

#dev.off()

6.1.3 Phylogenetic diversity

betadisper(beta_q1p$S, sample_metadata$location) %>% permutest(., pairwise=TRUE) 

Permutation test for homogeneity of multivariate dispersions
Permutation: free
Number of permutations: 999

Response: Distances
          Df Sum Sq Mean Sq      F N.Perm Pr(>F)  
Groups     5 0.2152 0.04304 2.1047    999  0.082 .
Residuals 86 1.7587 0.02045                       
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Pairwise comparisons:
(Observed p-value below diagonal, permuted p-value above diagonal)
               Aruba     Brazil  CaboVerde      Spain    Denmark Malaysia
Aruba                0.06500000 0.04800000 0.00200000 0.10900000    0.004
Brazil    0.06236518            0.85700000 0.33200000 0.66800000    0.679
CaboVerde 0.05279004 0.86912737            0.20700000 0.74600000    0.515
Spain     0.00060645 0.32391452 0.20693940            0.11000000    0.480
Denmark   0.09811567 0.63811338 0.73976492 0.10108967               0.300
Malaysia  0.00519311 0.67845753 0.52055144 0.50696000 0.31212601         
adonis2(beta_q1p$S ~ location, 
        data = sample_metadata %>% arrange(match(sample,labels(beta_q1p$S))), 
        permutations = 999) %>%
        broom::tidy() %>%
        tt()
tinytable_q1mfo9cgzd0fxnsj8wwf
term df SumOfSqs R2 statistic p.value
location 5 2.690791 0.2683828 6.309563 0.001
Residual 86 7.335152 0.7316172 NA NA
Total 91 10.025944 1.0000000 NA NA
pairwise.adonis(beta_q1p$S, sample_metadata$location, perm = 999)
                   pairs Df  SumsOfSqs    F.Model         R2 p.value p.adjusted sig
1        Aruba vs Brazil  1 0.64714035  6.1133065 0.16928127   0.001      0.015   .
2     Aruba vs CaboVerde  1 0.56229778  5.3504131 0.15575979   0.001      0.015   .
3       Aruba vs Denmark  1 1.03758657  9.5591181 0.24790811   0.001      0.015   .
4      Aruba vs Malaysia  1 0.48047063  5.0352414 0.14794199   0.001      0.015   .
5         Aruba vs Spain  1 0.60047832  6.7919744 0.18976250   0.001      0.015   .
6    Brazil vs CaboVerde  1 0.70961443  8.1023431 0.21837821   0.001      0.015   .
7      Brazil vs Denmark  1 0.32742702  3.5968626 0.11034383   0.014      0.210    
8     Brazil vs Malaysia  1 0.19570506  2.5119802 0.07971509   0.031      0.465    
9        Brazil vs Spain  1 0.06456161  0.9106384 0.03044530   0.435      1.000    
10  CaboVerde vs Denmark  1 1.07120011 11.9405505 0.29895808   0.001      0.015   .
11 CaboVerde vs Malaysia  1 0.50032095  6.5728056 0.19011490   0.002      0.030   .
12    CaboVerde vs Spain  1 0.82104026 11.9236999 0.29866220   0.001      0.015   .
13   Denmark vs Malaysia  1 0.54258755  6.8084757 0.19559821   0.001      0.015   .
14      Denmark vs Spain  1 0.24626102  3.3999417 0.10827860   0.017      0.255    
15     Malaysia vs Spain  1 0.26932060  4.5771932 0.14050299   0.001      0.015   .
#pdf("figures/beta_q1p_loca.pdf",width=9, height=5)
beta_q1n$S %>%
  vegan::metaMDS(., trymax = 500, k = 2, trace=0) %>%
  vegan::scores() %>%
  as_tibble(., rownames = "sample") %>%
  left_join(sample_metadata, by = "sample") %>%
  group_by(location) %>%
  mutate(x_cen = mean(NMDS1, na.rm = TRUE)) %>%
  mutate(y_cen = mean(NMDS2, na.rm = TRUE)) %>%
  ungroup() %>%
  ggplot(aes(x = NMDS1, y = NMDS2, color = location, fill = location)) +
    geom_point(size = 4) +
    #   stat_ellipse(aes(color = beta_q1n_nmds$Groups))+
    geom_segment(aes(x = x_cen, y = y_cen, xend = NMDS1, yend = NMDS2), alpha = 0.9, show.legend = FALSE) +
    scale_color_manual(values = location_colors)+
    theme_classic() +
    theme(
      axis.text.x = element_text(size = 12),
      axis.text.y = element_text(size = 12),
      axis.title = element_text(size = 20, face = "bold"),
      axis.text = element_text(face = "bold", size = 18),
      panel.background = element_blank(),
      axis.line = element_line(size = 0.5, linetype = "solid", colour = "black"),
      legend.text = element_text(size = 16),
      legend.title = element_text(size = 18),
      legend.position = "right", legend.box = "vertical"
    ) 

#dev.off()

6.1.4 Functional diversity

betadisper(beta_q1f$S, sample_metadata$location) %>% permutest(., pairwise=TRUE) 

Permutation test for homogeneity of multivariate dispersions
Permutation: free
Number of permutations: 999

Response: Distances
          Df Sum Sq  Mean Sq      F N.Perm Pr(>F)    
Groups     5 1.5016 0.300325 9.1732    999  0.001 ***
Residuals 86 2.8156 0.032739                         
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Pairwise comparisons:
(Observed p-value below diagonal, permuted p-value above diagonal)
               Aruba     Brazil  CaboVerde      Spain    Denmark Malaysia
Aruba                3.0000e-03 6.0600e-01 3.0000e-03 1.0000e-03    0.497
Brazil    8.6955e-04            2.0000e-03 6.9600e-01 2.4400e-01    0.003
CaboVerde 5.8936e-01 9.7727e-04            3.0000e-03 1.0000e-03    0.873
Spain     1.5125e-03 6.8358e-01 1.5804e-03            3.4000e-02    0.001
Denmark   1.7539e-04 2.3957e-01 9.2505e-05 3.3658e-02               0.001
Malaysia  4.9615e-01 1.0071e-03 8.8292e-01 1.5959e-03 7.4748e-05         
adonis2(beta_q1f$S ~ location, 
        data = sample_metadata %>% arrange(match(sample,labels(beta_q1f$S))), 
        permutations = 999) %>%
        broom::tidy() %>%
        tt()
tinytable_vmx95m7ty6n38djhpduh
term df SumOfSqs R2 statistic p.value
location 5 8.815599 0.5844961 24.19552 0.001
Residual 86 6.266793 0.4155039 NA NA
Total 91 15.082391 1.0000000 NA NA
pairwise.adonis(beta_q1f$S, sample_metadata$location, perm = 999)
                   pairs Df   SumsOfSqs    F.Model           R2 p.value p.adjusted sig
1        Aruba vs Brazil  1  3.77109355 40.6261637  0.575228238   0.001      0.015   .
2     Aruba vs CaboVerde  1  0.81599223  5.5801053  0.161367505   0.013      0.195    
3       Aruba vs Denmark  1  4.54372380 50.8948430  0.637022880   0.001      0.015   .
4      Aruba vs Malaysia  1  0.87672195  6.0397240  0.172367910   0.013      0.195    
5         Aruba vs Spain  1  3.49955224 36.9984243  0.560595570   0.001      0.015   .
6    Brazil vs CaboVerde  1  1.79573937 27.7696018  0.489163230   0.001      0.015   .
7      Brazil vs Denmark  1  0.27026868 35.0545409  0.547260825   0.001      0.015   .
8     Brazil vs Malaysia  1  1.69391038 26.6369446  0.478763612   0.001      0.015   .
9        Brazil vs Spain  1 -0.01258568 -0.9666596 -0.034482498   0.953      1.000    
10  CaboVerde vs Denmark  1  2.85656348 47.6207158  0.629731090   0.001      0.015   .
11 CaboVerde vs Malaysia  1 -0.02164970 -0.1836835 -0.006603446   0.978      1.000    
12    CaboVerde vs Spain  1  1.56212224 23.8545897  0.460028512   0.001      0.015   .
13   Denmark vs Malaysia  1  2.82610262 48.0023359  0.631590271   0.001      0.015   .
14      Denmark vs Spain  1  0.38029019 58.5480612  0.676480332   0.001      0.015   .
15     Malaysia vs Spain  1  1.47580731 22.9256154  0.450178466   0.001      0.015   .
#pdf("figures/beta_q1f_loca.pdf",width=9, height=5)
beta_q1f$S %>%
  vegan::metaMDS(., trymax = 500, k = 2, trace=0) %>%
  vegan::scores() %>%
  as_tibble(., rownames = "sample") %>%
  left_join(sample_metadata, by = "sample") %>%
  group_by(location) %>%
  mutate(x_cen = mean(NMDS1, na.rm = TRUE)) %>%
  mutate(y_cen = mean(NMDS2, na.rm = TRUE)) %>%
  ungroup() %>%
  ggplot(aes(x = NMDS1, y = NMDS2, color = location, fill = location)) +
    geom_point(size = 4) +
    #   stat_ellipse(aes(color = beta_q1n_nmds$Groups))+
    geom_segment(aes(x = x_cen, y = y_cen, xend = NMDS1, yend = NMDS2), alpha = 0.9, show.legend = FALSE) +
    scale_color_manual(values = location_colors)+
    theme_classic() +
    theme(
      axis.text.x = element_text(size = 12),
      axis.text.y = element_text(size = 12),
      axis.title = element_text(size = 20, face = "bold"),
      axis.text = element_text(face = "bold", size = 18),
      panel.background = element_blank(),
      axis.line = element_line(size = 0.5, linetype = "solid", colour = "black"),
      legend.text = element_text(size = 16),
      legend.title = element_text(size = 18),
      legend.position = "right", legend.box = "vertical"
    ) 

#dev.off()

6.2 Behaviour

6.2.1 Richness diversity

betadisper(beta_q0n$S, sample_metadata$origin) %>% permutest(., pairwise=TRUE) 

Permutation test for homogeneity of multivariate dispersions
Permutation: free
Number of permutations: 999

Response: Distances
          Df  Sum Sq   Mean Sq      F N.Perm Pr(>F)
Groups     1 0.00411 0.0041147 0.2442    999  0.613
Residuals 90 1.51634 0.0168483                     

Pairwise comparisons:
(Observed p-value below diagonal, permuted p-value above diagonal)
         Domestic Feral
Domestic          0.621
Feral     0.62238      
adonis2(beta_q0n$S ~ origin, 
        data = sample_metadata %>% arrange(match(sample,labels(beta_q0n$S))), 
        permutations = 999,
        strata = sample_metadata %>% arrange(match(sample,labels(beta_q1n$S))) %>% pull(location)) %>%
        broom::tidy() %>%
        tt()
tinytable_wjfsk1bqa1g86dnmdzrh
term df SumOfSqs R2 statistic p.value
origin 1 0.4407951 0.0163738 1.498173 0.178
Residual 90 26.4799587 0.9836262 NA NA
Total 91 26.9207538 1.0000000 NA NA
pairwise.adonis(beta_q0n$S, sample_metadata$origin, perm = 999)
              pairs Df SumsOfSqs  F.Model        R2 p.value p.adjusted sig
1 Domestic vs Feral  1 0.4407951 1.498173 0.0163738   0.081      0.081    
beta_q0n$S %>%
  vegan::metaMDS(., trymax = 500, k = 2, trace=0) %>%
  vegan::scores() %>%
  as_tibble(., rownames = "sample") %>%
  left_join(sample_metadata, by = join_by(sample == sample)) %>%
  group_by(origin) %>%
  mutate(x_cen = mean(NMDS1, na.rm = TRUE)) %>%
  mutate(y_cen = mean(NMDS2, na.rm = TRUE)) %>%
  ungroup() %>%
  ggplot(aes(x = NMDS1, y = NMDS2, color = origin, fill = origin)) +
    geom_point(size = 4) +
    scale_color_manual(values = origin_colors) +
    geom_segment(aes(x = x_cen, y = y_cen, xend = NMDS1, yend = NMDS2), alpha = 0.9, show.legend = FALSE) +
    theme_classic() +
    theme(
      axis.text.x = element_text(size = 12),
      axis.text.y = element_text(size = 12),
      axis.title = element_text(size = 20, face = "bold"),
      axis.text = element_text(face = "bold", size = 18),
      panel.background = element_blank(),
      axis.line = element_line(size = 0.5, linetype = "solid", colour = "black"),
      legend.text = element_text(size = 16),
      legend.title = element_text(size = 18),
      legend.position = "right", legend.box = "vertical"
    ) 

6.2.2 Neutral diversity

betadisper(beta_q1n$S, sample_metadata$origin) %>% permutest(., pairwise=TRUE) 

Permutation test for homogeneity of multivariate dispersions
Permutation: free
Number of permutations: 999

Response: Distances
          Df Sum Sq   Mean Sq      F N.Perm Pr(>F)
Groups     1 0.0080 0.0080048 0.5581    999  0.459
Residuals 90 1.2909 0.0143435                     

Pairwise comparisons:
(Observed p-value below diagonal, permuted p-value above diagonal)
         Domestic Feral
Domestic          0.463
Feral     0.45698      
adonis2(beta_q1n$S ~ origin, 
        data = sample_metadata %>% arrange(match(sample,labels(beta_q1n$S))), 
        permutations = 999,
        strata = sample_metadata %>% arrange(match(sample,labels(beta_q1n$S))) %>% pull(location)) %>%
        broom::tidy() %>%
        tt()
tinytable_m8sys3n5iaryj5ue10pr
term df SumOfSqs R2 statistic p.value
origin 1 0.386668 0.01498735 1.369385 0.315
Residual 90 25.412962 0.98501265 NA NA
Total 91 25.799630 1.00000000 NA NA
pairwise.adonis(beta_q1n$S, sample_metadata$origin, perm = 999)
              pairs Df SumsOfSqs  F.Model         R2 p.value p.adjusted sig
1 Domestic vs Feral  1  0.386668 1.369385 0.01498735    0.14       0.14    
pdf("figures/beta_q1n_behaviour.pdf",width=9, height=5)
beta_q1n$S %>%
  vegan::metaMDS(., trymax = 500, k = 2, trace=0) %>%
  vegan::scores() %>%
  as_tibble(., rownames = "sample") %>%
  left_join(sample_metadata, by = join_by(sample == sample)) %>%
  group_by(origin) %>%
  mutate(x_cen = mean(NMDS1, na.rm = TRUE)) %>%
  mutate(y_cen = mean(NMDS2, na.rm = TRUE)) %>%
  ungroup() %>%
  ggplot(aes(x = NMDS1, y = NMDS2, color = origin, fill = origin)) +
    geom_point(size = 4) +
    scale_color_manual(values = origin_colors) +
    geom_segment(aes(x = x_cen, y = y_cen, xend = NMDS1, yend = NMDS2), alpha = 0.9, show.legend = FALSE) +
    theme_classic() +
    theme(
      axis.text.x = element_text(size = 12),
      axis.text.y = element_text(size = 12),
      axis.title = element_text(size = 20, face = "bold"),
      axis.text = element_text(face = "bold", size = 18),
      panel.background = element_blank(),
      axis.line = element_line(size = 0.5, linetype = "solid", colour = "black"),
      legend.text = element_text(size = 16),
      legend.title = element_text(size = 18),
      legend.position = "right", legend.box = "vertical"
    ) 
dev.off()
RStudioGD 
        2 

6.2.3 Phylogenetic diversity

betadisper(beta_q1p$S, sample_metadata$origin) %>% permutest(., pairwise=TRUE) 

Permutation test for homogeneity of multivariate dispersions
Permutation: free
Number of permutations: 999

Response: Distances
          Df  Sum Sq  Mean Sq      F N.Perm Pr(>F)
Groups     1 0.01229 0.012295 0.5442    999  0.464
Residuals 90 2.03337 0.022593                     

Pairwise comparisons:
(Observed p-value below diagonal, permuted p-value above diagonal)
         Domestic Feral
Domestic          0.467
Feral     0.46262      
adonis2(beta_q1p$S ~ origin, 
        data = sample_metadata %>% arrange(match(sample,labels(beta_q1p$S))), 
        permutations = 999,
        strata = sample_metadata %>% arrange(match(sample,labels(beta_q1n$S))) %>% pull(location)) %>%
        broom::tidy() %>%
        tt()
tinytable_3uz62gakb3gwdppoveoq
term df SumOfSqs R2 statistic p.value
origin 1 0.2191073 0.02185403 2.010807 0.16
Residual 90 9.8068363 0.97814597 NA NA
Total 91 10.0259435 1.00000000 NA NA
pairwise.adonis(beta_q1p$S, sample_metadata$origin, perm = 999)
              pairs Df SumsOfSqs  F.Model         R2 p.value p.adjusted sig
1 Domestic vs Feral  1 0.2191073 2.010807 0.02185403   0.073      0.073    
beta_q1p$S %>%
  vegan::metaMDS(., trymax = 500, k = 2, trace=0) %>%
  vegan::scores() %>%
  as_tibble(., rownames = "sample") %>%
  left_join(sample_metadata, by = join_by(sample == sample)) %>%
  group_by(origin) %>%
  mutate(x_cen = mean(NMDS1, na.rm = TRUE)) %>%
  mutate(y_cen = mean(NMDS2, na.rm = TRUE)) %>%
  ungroup() %>%
  ggplot(aes(x = NMDS1, y = NMDS2, color = origin, fill = origin)) +
    geom_point(size = 4) +
    scale_color_manual(values = origin_colors) +
    geom_segment(aes(x = x_cen, y = y_cen, xend = NMDS1, yend = NMDS2), alpha = 0.9, show.legend = FALSE) +
    theme_classic() +
    theme(
      axis.text.x = element_text(size = 12),
      axis.text.y = element_text(size = 12),
      axis.title = element_text(size = 20, face = "bold"),
      axis.text = element_text(face = "bold", size = 18),
      panel.background = element_blank(),
      axis.line = element_line(size = 0.5, linetype = "solid", colour = "black"),
      legend.text = element_text(size = 16),
      legend.title = element_text(size = 18),
      legend.position = "right", legend.box = "vertical"
    ) 

6.2.4 Functional diversity

betadisper(beta_q1f$S, sample_metadata$origin) %>% permutest(., pairwise=TRUE) 

Permutation test for homogeneity of multivariate dispersions
Permutation: free
Number of permutations: 999

Response: Distances
          Df Sum Sq  Mean Sq      F N.Perm Pr(>F)  
Groups     1 0.2486 0.248604 3.3502    999  0.088 .
Residuals 90 6.6786 0.074206                       
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Pairwise comparisons:
(Observed p-value below diagonal, permuted p-value above diagonal)
         Domestic Feral
Domestic          0.091
Feral    0.070508      
adonis2(beta_q1f$S ~ origin, 
        data = sample_metadata %>% arrange(match(sample,labels(beta_q1f$S))), 
        permutations = 999,
        strata = sample_metadata %>% arrange(match(sample,labels(beta_q1n$S))) %>% pull(location)) %>%
        broom::tidy() %>%
        tt()
tinytable_1bqrpj7r7lh4h69mrtly
term df SumOfSqs R2 statistic p.value
origin 1 0.3816603 0.02530502 2.336579 0.18
Residual 90 14.7007312 0.97469498 NA NA
Total 91 15.0823914 1.00000000 NA NA
pairwise.adonis(beta_q1f$S, sample_metadata$origin, perm = 999)
              pairs Df SumsOfSqs  F.Model         R2 p.value p.adjusted sig
1 Domestic vs Feral  1 0.3816603 2.336579 0.02530502   0.127      0.127    
beta_q1f$S %>%
  vegan::metaMDS(., trymax = 500, k = 2, trace=0) %>%
  vegan::scores() %>%
  as_tibble(., rownames = "sample") %>%
  left_join(sample_metadata, by = join_by(sample == sample)) %>%
  group_by(origin) %>%
  mutate(x_cen = mean(NMDS1, na.rm = TRUE)) %>%
  mutate(y_cen = mean(NMDS2, na.rm = TRUE)) %>%
  ungroup() %>%
  ggplot(aes(x = NMDS1, y = NMDS2, color = origin, fill = origin)) +
    geom_point(size = 4) +
    scale_color_manual(values = origin_colors) +
    geom_segment(aes(x = x_cen, y = y_cen, xend = NMDS1, yend = NMDS2), alpha = 0.9, show.legend = FALSE) +
    theme_classic() +
    theme(
      axis.text.x = element_text(size = 12),
      axis.text.y = element_text(size = 12),
      axis.title = element_text(size = 20, face = "bold"),
      axis.text = element_text(face = "bold", size = 18),
      panel.background = element_blank(),
      axis.line = element_line(size = 0.5, linetype = "solid", colour = "black"),
      legend.text = element_text(size = 16),
      legend.title = element_text(size = 18),
      legend.position = "right", legend.box = "vertical"
    )